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ABSTRACT 


The  periodic  shedding  of  vortices  from  bluff  bodies  was  first  recognized  in  the 
late  1800’s.  Currently,  there  is  great  interest  concerning  the  effect  of  vortex  shedding 
on  structures  and  on  vehicle  stability.  In  the  design  of  bluff  structures  which  will  be 
exposed  to  a  flow,  knowledge  of  the  shedding  frequency  and  the  amplitude  of  the 
aerodynamic  forces  is  critical.  The  ability  to  computationally  predict  parameters 
associated  with  periodic  vortex  shedding  is  thus  a  valuable  tool. 

In  this  study,  the  periodic  shedding  of  vortices  from  several  bluff  body  geome¬ 
tries  is  predicted.  The  study  is  conducted  with  a  two-dimensional  finite-difference 
code  employed  on  various  grid  sizes.  The  effects  of  the  grid  size  and  time  step  on 
the  accuracy  of  the  solution  are  addressed.  Strouhal  numbers  and  aerodynamic 
force  coefficients  are  computed  for  all  of  the  bodies  considered  and  compared  with 
previous  experimental  results. 

Results  indicate  that  the  finite-difference  code  is  capable  of  predicting  periodic 
vortex  shedding  for  all  of  the  geometries  tested.  Refinement  of  the  finite-difference 
grid  was  found  to  give  little  improvement  in  the  prediction;  however,  the  choice  of 
time  step  size  was  shown  to  be  critical.  Predictions  of  Strouhal  numbers  were  gen¬ 
erally  accurate,  and  the  calculated  aerodynamic  forces  generally  exhibited  behavior 
consistent  with  previous  studies.  •  \  » 
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ABSTRACT 

The  periodic  shedding  of  vortices  from  bluff  bodies  was  first  recognized  in  the 
late  1800’s.  Currently,  there  is  great  interest  concerning  the  effect  of  vortex  shedding 
on  structures  and  on  vehicle  stability.  In  the  design  of  bluff  structures  which  will  be 
exposed  to  a  flow,  knowledge  of  the  shedding  frequency  and  the  amplitude  of  the 
aerodynamic  forces  is  critical.  The  ability  to  computationally  predict  parameters 
associated  with  periodic  vortex  shedding  is  thus  a  valuable  tool. 

In  this  study,  the  periodic  shedding  of  vortices  from  several  bluff  body  geome¬ 
tries  is  predicted.  The  study  is  conducted  with  a  two-dimensional  finite-difference 
code  employed  on  various  grid  sizes.  The  effects  of  the  grid  size  and  time  step  on 
the  accuracy  of  the  solution  are  addressed.  Strouhal  numbers  and  aerodynamic 
force  coefficients  are  computed  for  all  of  the  bodies  considered  and  compared  with 
previous  experimental  results. 

Results  indicate  that  the  finite-difference  code  is  capable  of  predicting  periodic 
vortex  shedding  for  all  of  the  geometries  tested.  Refinement  of  the  finite-difference 
grid  was  found  to  give  little  improvement  in  the  prediction;  however,  the  choice  of 
time  step  size  was  shown  to  be  critical.  Predictions  of  Strouhal  numbers  were  gen¬ 
erally  accurate,  and  the  calculated  aerodynamic  forces  generally  exhibited  behavior 
consistent  with  previous  studies. 
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NOMENCLATURE 


constants 


ag,  O/v,  as,  a\v  finite-difference  coefficients 
A  frontal  area 


source  term 


length  of  rectangular  cylinder 

nondimensional  drag  coefficient 

nondimensional  lift  coefficient 

pressure  coefficient 

base  pressure  coefficient 

characteristic  length 

wake  width 

drag  force 

internal  energy 

function  dependent  on  A 

height  of  plate  above  the  ground  plane 

kinetic  energy  of  turbulence 

plate  or  cylinder  height 


lift  force 


number  of  cells  in  a  horizontal  wall 

shedding  frequency 

number  of  cells  in  a  vertical  wall 

static  pressure 

freestream  static  pressure 

shedding  period 

dynamic  pressure 

added  heat 

specific  gas  constant 

Reynolds  number 

Strouhal  number 

function  dependent  on  Re 

time 

time  step 
temperature 
x  component  of  velocity 
relative  velocity  of  vortices 
separation  velocity 
freestream  velocity 
y  component  of  velocity 
velocity  vector 
z  component  of  velocity 
work  done  by  the  system 
x  direction  body  force 


y  direction  body  force 
z  direction  body  force 
relaxation  factor 

aspect  ratio  of  rectangular  cylinder 

dissipation  rate  of  turbulence  energy 

ratio  of  Us  to 

ratio  of  h  to  / 

dynamic  viscosity 

density 

shear  stress 

general  variable 


CHAPTER  I 
INTRODUCTION 


Formal  recognition  of  the  periodic  structures  existing  in  the  wakes  of  bluff 
bodies  is  credited  to  Strouhal.  In  his  now  famous  work  of  1878,  Strouhal  determined 
that  the  shedding  of  vortices  behind  a  body  occurred  periodically,  and  he  also  gave 
a  crude  approximation  for  the  frequency  of  the  vortex  shedding.  Soon  after  in  1879, 
Rayleigh  discovered  a  relationship  between  the  shedding  frequency  and  the  flow 
Reynolds  number.  However,  the  most  famous  contributions  to  the  field  of  periodic 
structures  in  bluff  body  wakes  is  credited  to  von  Karman.  Von  Karman  formulated 
a  theory  on  the  flow  pattern  that  occurs  in  the  wakes  of  bluff  bodies,  that  pattern 
being  the  familiar  Karman  vortex  street  (Figure  1.1).  His  theory  attempted  to 
predict  the  geometry  of  the  flow  pattern  in  the  wake  behind  a  body  as  well  as  other 
vital  parameters,  such  as  the  velocity  of  the  vortices  and  the  drag  on  the  body. 
Though  his  theory  is  not  complete  in  that  it  concentrates  only  on  the  wake  region 
and  subsequently  ignores  the  exterior  potential  flow,  it  is  nonetheless  a  significant 
contribution  to  the  study  of  coherent  structures  in  bluff  body  wakes. 


Importance  of  Vortex  Shedding 


These  pioneering  efforts  in  the  study  of  periodic  structures  in  wakes  have  many 
valuable  applications.  Though  the  study  of  the  vortex  street  is  essentially  a  problem 
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in  fluid  dynamics,  its  primary  applications  are  related  to  structural  mechanics.  Since 
the  vortex  street  is  periodic  in  nature,  the  aerodynamic  loads  on  a  body  which 
is  shedding  vortices  are  also  periodic.  The  frequency  of  the  vortex  shedding  is  a 
critical  parameter,  for  if  the  shedding  frequency  approaches  a  natural  frequency 
of  the  structure,  potentially  damaging  large  amplitude  oscillations  can  occur.  A 
vivid  demonstration  of  the  damaging  effects  of  vortex  shedding  was  provided  by 


the  destruction  of  the  Tacoma  Narrows  Bridge  in  1940  (Polodny  and  Scalzi  1976). 
In  a  matter  of  hours,  this  bridge  was  destroyed  by  the  oscillations  induced  by  the 
periodic  shedding  of  vortices  from  the  bridge’s  surface. 

In  line  with  this  example,  the  study  of  vortex  shedding  is  important  in  the 
design  of  structures  that  will  be  subjected  to  wind  loadings.  Also,  there  is  an  analog 
in  marine  structures,  where  vortex  shedding  is  a  matter  of  considerable  interest. 
Another  field  of  recent  interest  is  the  effect  of  vortex  shedding  on  the  stability  of 
motor  vehicles  with  bluff  tail  sections.  In  general,  any  bluff  body  subjected  to  a 
fluid  flow  can  potentially  shed  vortices;  knowledge  of  the  frequency  and  amplitude 
of  the  associated  loading  can  be  vital  to  stability  and  structural  integrity. 


The  ability  to  predict  various  parameters  associated  with  vortex  shedding  can 
be  of  great  value.  For  structures  which  will  be  subjected  to  fluid  flow,  the  possibility 
of  damage  due  to  vortex -induced  loading  must  be  considered.  Though  one  could 
attempt  to  measure  the  vortex  shedding  characteristics  experimentally,  the  ability 
to  do  so  computationally  has  many  advantages. 

The  advantages  of  computational  studies  are  often  based  on  monetary  consid¬ 
erations.  One  fundamental  advantage  of  computational  methods  is  the  ability  to 


construct  and  test  many  configurations  relatively  quickly  and  inexpensively.  The 
use  of  fine  grids  in  calculations  allows  one  to  achieve  a  highly  detailed  picture  of 
the  flow  field,  and  greater  detail  can  be  achieved  by  simply  refining  the  grid.  In  an 


experimental  study,  determining  complex  flow  fields  may  require  time  consuming 
testing,  which  inevitably  leads  to  great  expense.  In  many  cases,  computer  codes 
have  the  ability  to  predict  flow  in  regions  of  a  flow  field  that  are  inaccessible  by  ex¬ 
perimental  methods. 

This  study  treats  the  problem  of  vortex  shedding  with  a  two-dimensional  finite- 
difference  code.  A  three-dimensional  code  is  not  used  for  several  reasons.  First, 
there  was  no  three-dimensional  code  available  for  this  study.  Also,  the  addition  of 
another  dimension  to  the  calculation  domain  would  greatly  increase  computer  run 
time.  Furthermore,  since  the  primary  objective  of  the  study  is  to  simply  demon¬ 
strate  the  ability  to  numerically  predict  vortex  shedding,  the  added  complexity  of 
a  third  dimension  is  unnecessary. 

Objectives  of  the  Study 

The  primary  objective  of  this  study  is  to  demonstrate  the  ability  of  a  finite- 
difference  code  to  predict  periodic  vortex  shedding  from  various  geometries.  To 
demonstrate  the  abilities  of  the  code,  the  phenomena  of  vortex  shedding  is  consid¬ 
ered  for  three  basic  flow  schemes.  The  three  basic  schemes  consist  of  the  following: 
a  flat  plate  in  a  uniform  flow  which  is  oriented  normal  to  the  plate  (Figure  1.2(a)), 
a  flat  plate  near  a  ground  plane  in  a  uniform  flow  which  is  oriented  normal  to  the 
plate  (Figure  1.2(b)),  and  a  rectangular  cylinder  in  a  uniform  flow  which  is  oriented 
normal  to  the  cylinder  (Figure  1.2(c)).  The  dependence  of  the  results  on  the  choice 
of  finite-difference  grid  size,  turbulence,  and  time  step  size  are  also  examined. 
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Data  for  the  comparison  of  computational  and  experimental  results  comes  from 
various  sources.  For  the  shedding  from  a  flat  plate  in  a  open  region,  the  parameters 
used  for  the  comparison  are  shedding  frequency,  drag  and  lift  coefficients,  and  wake 
pressure  profiles.  Experimental  data  for  comparison  comes  primarily  from  work 
by  Roshko  (1954).  For  the  case  of  a  flat  plate  near  a  ground  plane,  the  primary 
interest  is  the  prediction  of  shedding  frequency.  Experimental  results  are  provided 
by  the  studies  of  Strickland  et  al.  (1980)  and  Matty  (1979).  Finally,  the  study  of 
vortex  shedding  from  a  rectangular  body  concerns  the  prediction  of  the  shedding 
frequency  and  the  drag  coefficient.  Data  for  comparison  is  extracted  from  work  by 
Okajima  (1982)  and  Davis  and  Moore  (1982). 


(a)  Flat  Plate  in  a  Uniform  Flow. 


(b)  Flat  Plate  Near  a  Ground  Plane 
in  a  Uniform  Flow. 


(c)  Rectangular  Cylinder  in  a 
Uniform  Flow. 


FIGURE  1.2:  Basic  Flow  Configurations. 
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CHAPTER  II 


THEORY  OF  VORTEX  SHEDDING 

The  subsequent  discussion  covers  the  aspects  of  vortex  shedding  which  are 
relevant  to  the  current  study.  The  basic  mechanism  of  periodic  vortex  shedding  is 
discussed  along  with  a  review  of  experimental  results  for  each  of  the  three  geometries 
studied.  The  experimental  results  presented  provide  a  basis  for  comparison  with 
the  numerical  p  idictions. 

Basic  Parameters 

The  study  of  vortex  shedding  is  simplified  by  the  introduction  of  a  number  of 
dimensionless  parameters.  Several  parameters  are  defined  to  nondimensionalize  the 
aerodynamic  forces,  pressure,  and  vortex  shedding  frequency. 

The  fundamental  aerodynamic  forces  known  as  drag  and  lift  are  typically  given 
in  a  nondimensional  form.  The  nondimensional  forms  of  drag  and  lift  are  respec¬ 
tively  the  drag  coefficient,  Cd,  and  the  lift  coefficient,  Cl-  Definitions  of  Cd  and 
Cl  are  given  respectively  as 


Cd  = 


D 


§00-4 


(2.1) 

(2.2) 
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=  frontal  area 
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Z)  =  drag  force 

L  =  lift  force 

<7oo  =  dynamic  pressure. 

The  dynamic  pressure  is  defined  as 


9°o  —  2^°°  (2-3) 

where 

U oo  =  freest  ream  velocity 

p  =  density. 

Pressures  are  nondimensionalized  by  introducing  the  pressure  coefficient,  Cp, 
defined  as 


cp  = 


P  -  Poo 


(2.4) 
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where 

p  =  static  pressure 

Poo  =  freestream  static  pressure. 

Note  that  Cp  >  0  indicates  that  the  pressure  is  greater  than  freestream  pressure, 
and  Cp  <  0  indicates  that  the  pressure  is  less  than  freestream  pressure. 

Another  nondimensional  quantity  often  associated  with  flow  problems  is  the 
Reynolds  number,  Re.  The  Reynolds  number  expresses  the  ratio  of  inertial  to 
viscous  forces  and  is  defined  as 


2 


i 


where 


d  —  characteristic  length 

/x  =  dynamic  viscosity. 

The  final  parameter  to  introduce  is  the  Strouhal  number,  St.  The  Strouhal 
number  is  a  nondimensional  frequency,  and  it  is  the  primary  quantity  in  the  study 
of  vortex  shedding.  The  Strouhal  number  is  defined  as 


St  = 


nd 

uZ 


where  n  is  the  shedding  frequency. 


(2.6) 


Vortex  Formation  from  Bluff  Bodies 

A  bluff  body  in  a  uniform  flow  can  assume  two  different  flow  patterns  in  its 
wake.  At  very  low  Reynolds  numbers  (Re  <  60)  a  stationary,  symmetric  vortical 
pattern  is  formed.  This  symmetric  pattern  is  a  seemingly  trivial  solution,  because 
it  can  only  be  achieved  under  controlled  situations.  The  solution  of  interest  is  the 
formation  of  a  vortex  street  in  the  wake  of  the  body. 

The  mechanism  of  wake  formation  is  described  similarly  by  Younis  (1988)  and 
Roshko.  At  the  sharp  leading  edges  of  the  body  the  flow  separates,  creating  free 
shear  layers.  These  shear  layers  continue  downstream  for  a  short  distance  before 
rolling  up  into  vortices.  This  flow  pattern  is  inherently  unstable,  and  thus  oscilla¬ 
tions  in  the  flow  field  are  initiated.  Eventually,  the  vortices  from  the  top  and  bot¬ 
tom  surfaces  will  alternately  be  carried  downstream  forming  the  familiar  pattern  of 
the  Karman  vortex  street. 


Roshko  derived  certain  qualitative  relationships  pertaining  to  bluff  body  wakes. 
He  noted  that  all  bluff  body  wakes  are  similar  in  structure,  with  the  bluffness  of  a 
body  being  related  to  the  dimension  of  the  wake  width  with  respect  to  the  body.  In 
essence,  a  bluffer  body  diverges  the  flow  more  extensively  and  creates  a  larger  wake. 
He  also  noted  an  inverse  relationship  between  the  shedding  frequency  and  the  wake 
width,  indicating  that  bluffer  bodies  have  lower  shedding  frequencies.  Finally,  he 
noted  that  a  decrease  in  wake  width  is  accompanied  by  an  increased  drag. 

Vortex  Shedding  from  a  Flat  Plate 

The  shedding  of  vortices  from  bluff  bodies  was  studied  extensively  by  Roshko. 
Among  the  bodies  Roshko  studied  was  the  flat  plate.  Of  primary  interest  was  the 
measurement  of  shedding  frequency  from  the  plate  in  terms  of  the  Strouhal  number. 
He  determined  that  the  Strouhal  number  for  a  flat  plate  remained  nearly  constant 
over  a  wide  range  of  Reynolds  numbers.  The  value  of  the  Strouhal  number  for  a 
flat  plate  fell  in  the  range  of  0.130  <  St  <  0.140.  Experimental  Strouhal  numbers 
determined  by  Roshko  are  provided  as  Figure  2.1. 

Also  of  interest  was  the  prediction  of  aerodynamic  forces  on  the  plate,  especially 
drag.  Roshko  devised  a  method  by  which  he  could  determine  the  drag  on  a  flat 
plate  from  one  experimental  measurement,  the  base  pressure  coefficient.  He  joined 
the  work  of  Kirchhoff  and  von  Karman  in  what  was  termed  the  notched-hodograph 
theory  (Roshko  1953).  The  geometry  employed  in  the  notched-hodograph  theory  is 
presented  in  Figure  2.2,  with  the  free  streamlines  representing  the  separated  shear 
layers. 

The  velocity  Us  is  the  separation  velocity  and  is  a  multiple  of  the  freestream 


velocity  by  the  relation 
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lr„  =  ict/Q 


Since  the  velocity  U„  exists  on  the  downwind  side  of  the  plate,  the  base  pressure 


coefficient,  Cps,  is  simply 


Cps  =  1  -  k2  . 


To  determine  k,  the  relationship  is  transformed  to  read 


k  —  \/l  —  C  T 


The  constant  k  is  related  to  the  wake  width  by  the  curve  given  as  Figure  2.3. 
Another  relation  was  formed  between  the  relative  velocity  of  the  vortices,  urti,  and 
the  freestream  velocity  as 


UTel  1 


Finally,  the  drag  can  be  calculated  by  the  Karman  drag  formula 


(2.10) 


CD  =  j  5.65^  -2.25^ 


(2.11) 


The  data  provided  by  Roshko  for  Cp ,  yields  drag  coefficients  in  the  neighborhood 
of  1.85. 

Vortex  Sheddine  from  a  Flat  Plate 


Near  a  uround  Plane 

The  study  of  vortex  shedding  from  a  body  situated  near  a  ground  plane  seems 
to  have  been  exclusively  addressed  by  the  related  studies  of  Matty  and  Strickland 
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et  al.  Their  study  was  motivated  by  a  concern  about  the  effects  of  vortex  shedding 
from  heliostats.  Though  the  actual  geometry  that  they  were  concerned  with  was 
not  a  flat  plate,  they  addressed  their  concern  by  performing  wind  tunnel  tests  on 
a  flat  plate  near  a  ground  plane.  The  basic  configuration  involved  in  the  study  is 
shown  in  Figure  2.4. 

Their  study  intended  to  determine  the  Strouhal  number  as  a  function  of  the 
Reynolds  number  and  the  height  of  the  plate  above  the  ground  plane.  Shedding 
frequencies  were  determined  for  the  three-dimensional  flat  plate,  and  the  data  they 
obtained  are  presented  in  Figure  2.5.  From  these  data,  an  equation  was  derived 
that  approximates  the  Strouhal  number  in  terms  of  the  dimensionless  parameters 
Re  and  A,  where 

A  =  h/l  .  (2.12) 

The  equation  derived  from  the  study  took  the  following  form: 

St(Re,  A)  =  gil(\)St00(Re)  .  (2.13) 

Their  derivation  gave  the  Strouhal  number  as  the  product  of  two  indepen¬ 
dent  functions.  The  function  gh{ A)  describes  the  variation  in  Strouhal  number  for 
changes  in  A.  The  function  g^(A)  is  given  as 

gh  =  1  +  a1e_ajA  sin[a3(a4A  -  1)]  (2.14) 


where  are  constants  given  as 


a2  —  4.6 

a3  =  18.0° 

a4  =  20.0. 

A  graphical  representation  of  the  function  gk{ A)  is  provided  as  Figure  2.6. 

The  function  St<x(Re )  is  representative  of  the  variation  in  the  Strouhal  number 
with  Reynolds  number  for  a  flat  plate  removed  from  the  ground  plane.  Analytically, 
the  function  Strxi(Re )  is  given  by  the  simple  expression 


Stoo  =  a5  +  a6/ Re  (2.15) 

where 

a5  =  0.121 

a6  =  407. 

This  function  is  shown  graphically  as  Figure  2.7.  It  is  important  to  note  that  this 
function  was  derived  for  the  specific  Reynolds  number  range  of  1.4  x  104  <  Re  < 
1.5  x  10s. 

There  are  a  few  important  features  to  note  about  the  function  St(Re,  A).  In¬ 
creasing  Reynolds  number  is  accompanied  by  a  decrease  in  the  Strouhal  number  as 
indicated  by  St00(Re).  With  respect  to  the  function  for  0  <  A  <  0.20  the  ten¬ 

dency  is  for  the  Strouhal  number  to  increase  with  increasing  A.  For  0.20  <  A  <  0.50, 
the  Strouhal  number  decreases  for  increasing  A.  Finally,  for  A  >  0.50,  the  Strouhal 
number  appears  to  be  virtually  unaffected  by  the  presence  of  the  ground  plane. 
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automobiles,  buildings,  and  many  other  bodies  of  practical  interest.  These  bodies 
also  exhibit  interesting  behavior  with  relation  to  their  characteristics  while  shedding 
vortices. 


The  more  curious  aspects  of  vortex  shedding  from  rectangular  cylinders  are  re¬ 
lated  to  the  Strouhal  numbers  of  such  bodies.  The  basic  geometry  of  the  rectangu¬ 
lar  cylinder  is  provided  as  Figure  2.8.  The  ratio  of  c/l,  which  is  a  critical  parame¬ 
ter  in  determining  the  characteristics  of  the  flow  about  a  given  cylinder,  is  denoted 
as  (3.  Okajima  provided  plots  of  the  variation  in  Strouhal  number  with  Reynolds 
number  for  four  values  of  /3  ((3  =  1,  2,  3,  4).  These  plots  (Figures  2.9-2.12)  demon¬ 
strate  the  interesting  behavior  apparent  in  the  shedding  frequency. 

For  the  case  (3  =  1  (Figure  2.9),  there  is  some  variation  in  the  Strouhal  number 
for  increasing  Reynolds  number,  but  the  variation  is  gradual.  However,  for  the  cases 
of  (3  =  2  (Figure  2.10)  and  /3  =  3  (Figure  2.11),  the  behavior  is  quite  different.  For 
(3  =  2,  there  is  a  discontinuity  in  the  plot  at  Re  %  500  and  a  similar  discontinuity 
appears  in  the  plot  for  f3  =  3  at  Re  %  800.  This  discontinuity  appears  because  there 
are  actually  two  modes  of  vortex  shedding  which  occur  from  rectangular  cylinders 
in  the  range  of  2  <  (3  <  3  (Okajima). 

The  first  mode,  which  occurs  at  low  Reynolds  numbers,  is  characterized  by  a 
reattachment  of  the  flow  which  separates  at  the  leading  edge.  The  reattached  flow- 
then  separates  at  the  trailing  edge  and  forms  the  vortex  street.  As  the  Reynolds 
number  is  increased,  the  flow  loses  the  ability  to  reattach  to  the  surface  of  the 
cylinder.  Finally,  when  the  Reynolds  number  reaches  a  critical  value,  the  separated 
regions  suddenly  detach  from  the  surfaces  creating  a  drastic  increase  in  the  wake 
width.  As  demonstrated  by  Roshko,  an  increase  in  the  wake  width  is  accompanied 


by  a  decrease  in  the  shedding  frequency.  Therefore,  this  detachment  characterizes 
the  second  mode  of  the  vortex  shedding  and  is  the  physical  process  responsible 
for  the  discontinuities  which  occur  in  plots  for  Strouhal  number  versus  Reynolds 
number. 

With  respect  to  the  variation  in  aerodynamic  forces  on  rectangular  cylinders, 
some  general  trends  have  been  recognized.  A  study  by  Davis  and  Moore  indicated 
an  increase  in  lift  and  drag  for  an  increase  in  Reynolds  number.  Laneville  and  Yong 
(1983)  were  concerned  with  the  variation  in  the  drag  for  changing  0.  Their  data 
indicate  that  the  drag  increases  steadily  with  increasing  0  for  0  <  0.8.  However,  for 
0  >  0.8  the  trend  suddenly  changes  and  increasing  0  is  accompanied  by  a  decrease 
in  the  drag  coefficient. 
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FIGURE  2.1:  Strouhal  Number  vs.  Reynolds  Number 
for  a  Flat  Plate  (Roshko  1954). 
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FIGURE  2.2:  Geometry  for  the  Notched-Hodograph  Theory. 
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FIGURE  2.5:  Strouhal  Number  vs.  Reynolds  Number  for 
a  Flat  Plate  Near  a  Ground  Plane. 


FIGURE  2.6  The  Function  gh( A)  (Matty  1979) 


FIGURE  2.11:  Strouhal  Number  vs.  Reynolds  Number 

for  a  Cylinder  with  0  =  3  (Okajima  1982) 


FIGURE  2.12:  Strouhal  Number  vs.  Reynolds  Number 

for  a  Cylinder  with  0  =  4  (Okajima  1982). 
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DESCRIPTION  OF  COMPUTATIONAL  ANALYSIS 
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The  field  of  computational  fluid  dynamics  is  currently  a  field  of  considerable 
interest.  The  basic  goal  is  to  predict  the  distributions  of  velocity,  pressure,  tem¬ 
perature  and  other  relevant  variables  throughout  a  given  flow  field.  Advantages  of 
computational  methods  over  experimental  methods  include  a  greatly  reduced  cost 
and  the  ability  to  survey  areas  of  a  flow  field  that  may  be  otherwise  inaccessible. 
Furthermore,  numerical  methods  allow  the  solution  of  problems  which  are  impos¬ 
sible  to  solve  using  classical  methods.  Computational  methods,  of  course,  are  not 
without  disadvantages.  In  some  cases,  numerical  solutions  may  diverge  and  yield 
no  useful  results.  Also,  for  more  complex  flows,  especially  turbulent  flows,  empiri¬ 
cal  models  are  employed  that  may  not  suitably  represent  the  physics  of  a  situation. 
This  may  lead  one  to  question  the  validity  of  computational  results.  For  this  reason, 
a  combination  of  experimental  and  computational  investigations  will  likely  provide 
the  most  consistent  and  useable  results.  In  the  following  discussion,  the  computa¬ 
tional  tools  employed  in  this  research  are  reviewed. 


Discretization  of  Differential  Equations 


For  a  homogeneous  fluid,  a  total  description  of  the  flow  field  consists  of  the 
distribution  of  six  quantities  throughout  the  flow  domain.  In  Cartesian  coordinates, 
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the  six  quantities  which  arise  in  these  equations  are  the  three  velocity  components, 
u,  v,  and  w,  pressure  p,  density  p  and  temperature  T .  The  required  six  equations 
for  the  solution  of  any  flow  problem  are  fundamental.  Three  equations  are  provided 
by  conservation  of  momentum  for  each  of  the  three  Cartesian  coordinate  directions 
( x,y,z ).  These  equations  are  also  commonly  known  as  the  Navier-Stokes  equations 
and  are  given  for  the  x,  y ,  and  z  directions  respectively  as 
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t  —  time 

X  —  i-direction  body  force 

Y  =  y-direction  body  force 

Z  =  2-direction  body  force 

fj.  =  dynamic  viscosity 
and  v  is  the  velocity  vector 


v  =  ui  +  uj  +  wk  .  (3-4) 

A  fourth  equation  is  provided  by  conservation  of  mass,  or  the  continuity  equation, 
given  as 


dp  d(pu )  8(pv)  d(pw) 

dt  dx  dy  dz 


(3.5) 


Conservation  of  energy  provides  the  fifth  equation  which  given  in  a  basic  form  is 


SQ  =  dET  +  6W  (3.6) 

where 

8Q  =  quantity  of  heat  added 
dEr  —  increase  in  internal  energy 
8W  =  work  done  by  the  system. 

The  set  of  equations  is  closed  by  an  equation  of  state,  commonly  the  well  known 
ideal  gas  law 

p  =  pRT  (3.7) 
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where  R  is  the  specific  gas  constant.  Of  course,  the  use  of  the  ideal  gas  equation 
assumes  that  the  fluid  is  a  gas.  If  this  were  not  the  case,  an  appropriate  equation 
of  state  would  have  to  be  employed. 

Theoretically,  given  six  equations  for  six  unknowns  the  solution  to  a  given  prob¬ 
lem  should  be  within  reach.  However,  with  the  exception  of  very  simple  problems, 
this  set  of  equations  cannot  be  solved  analytically.  It  then  becomes  the  task  of  the 
computational  fluid  dynamics  algorithm  to  provide  an  approximate  solution  that 
adequately  describes  the  flow  field. 

To  obtain  a  solution  to  a  given  flow  field  using  a  finite-difference  scheme,  the 
differential  equations  that  describe  the  physics  of  the  flow  must  first  be  discretized. 
To  discretize  the  differential  equations,  one  assumes  that  the  flow  field  is  composed 
of  a  discrete  number  of  control  volumes,  and  the  flow  variables  stored  at  each 
control  volume  are  prevalent  for  that  region  in  space.  The  discretization  equations 
are  algebraic  equations  that  relate  the  quantities  at  a  grid  location  to  the  quantities 
at  the  neighboring  grid  locations.  The  equations  that  are  derived  to  represent 
conservation  of  mass,  momentum,  and  energy  all  follow  the  same  basic  form,  which  is 

ap<pp  =  0L£(f>E  +  aw<}>w  +  aN4>N  +  as(f>s  +  b  (3-8) 

where 

a  =  finite-difference  coefficient 

b  =  source  term 

< I >  =  general  variable  (u,  v,  T,  etc.). 

The  subscripts  identify  which  point  in  the  field  the  a  or  <j>  refers  to,  as  shown  in 
Figure  3.1.  The  source  term  b  is  available  to  account  for  any  terms  which  do  not 
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conform  readily  to  this  format  upon  discretization  of  the  differential  equation.  The 
most  obvious  example  of  a  term  that  is  contained  in  the  source  term  b  is  the  pressure 
gradient  term  in  the  momentum  equations. 

Among  the  most  important  aspects  in  the  derivation  of  the  finite-difference 
coefficients  is  the  differencing  scheme.  The  differencing  scheme  provides  information 
on  the  variation  of  flow  quantities  between  grid  locations.  The  scheme  employed  in 
tills  research  is  known  as  the  upwind  scheme.  In  this  scheme,  the  value  of  a  quantity 
at  a  given  point  is  assumed  to  be  affected  only  by  the  cells  which  are  upstream  of 
the  cell,  with  no  influence  of  the  downstream  cells  being  felt.  This  scheme  works 
well  for  processes  governed  primarily  by  convection,  but  it  is  not  suitable  for  flows 
which  are  predominantly  processes  of  diffusion.  The  solution  of  these  discretization 
equations  throughout  the  flow  domain  is  the  task  of  the  finite-difference  algorithm. 
For  a  detailed  description  of  the  finite-difference  coefficients  and  other  differencing 
schemes  refer  to  Patankar  (1980). 

The  FORDC-2  Finite-Difference  Code 

The  basis  for  the  code  used  in  this  research  is  the  FORDC-2  (Flow  Over 
Radiator  and  Condenser  —  2  Dimensional)  code  developed  by  Carroll,  Maxwell  and 
Sun  (1985).  This  code  originally  treated  a  flow  as  two-dimensional,  incompressible, 
steady  and  elliptic.  Modifications  have  been  made  to  the  code  to  allow  the  solution 
of  unsteady  problems  and  to  speed  the  solution  of  the  discretization  equations. 

The  algorithm  employed  in  the  FORDC-2  code  is  the  SIMPLE  (Semi-  Implicit 
Method  for  Pressure  Linked  Equations)  algorithm  described  by  Patankar  (1980). 
A  brief  outline  of  this  iterative  method  follows.  The  iterative  process  is  started  by 
assuming  a  pressure  field.  With  the  pressure  field  presumed  to  be  known,  one  can 
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solve  the  momentum  equations  to  obtain  the  velocities.  Next,  the  pressure  correc¬ 
tion  equation,  a  form  of  the  continuity  equation  employed  to  insure  conservation  of 
mass,  is  solved.  Following  the  solution  of  the  pressure  correction  equation,  the  ve¬ 
locities  and  pressures  are  adjusted  to  insure  satisfaction  of  continuity.  The  solution 
then  proceeds  to  deal  with  any  other  quantities  of  interest  (i.e.  temperature).  The 
updated  pressure  is  then  used  to  continue  the  iterative  process.  Typically,  the  iter¬ 
ative  process  is  continued  until  the  flow  quantities  remain  virtually  unchanged  be¬ 
tween  iterations  indicating  a  converged  solution.  However,  in  the  transient  case  the 
solution  should  be  valid  at  each  step  in  time.  Flow  quantities  may  well  be  changing 
considerably  between  successive  iterations,  or  time  steps,  in  transient  problems. 

The  FORDC-2  code  was  modified  to  speed  the  solution  by  altering  the  method 
by  which  the  discretization  equations  are  solved.  Previously,  all  of  the  discretization 
equations  except  pressure  correction  were  solved  by  the  ADI  (Alternating  Direction 
Iteration)  procedure.  The  program  has  been  altered  to  allow  the  solution  of  all  of  the 
discretization  equations  by  Stone’s  method.  Stone’s  method  is  a  field  iterative  solver 
that  is  claimed  to  be  approximately  three  times  as  fast  as  other  available  meth¬ 
ods.  Stone’s  method  requires  the  selection  of  a  relaxation  factor  which  was  set  to 
a  =  .93.  For  details  on  the  formulation  of  Stone’s  method,  refer  to  either  Stone 
(1968)  or  Carroll  (1980). 

Finite-Difference  Grid 

The  accuracy  of  the  solution  of  the  finite-difference  equations  is  governed  by 
the  approximation  involved  in  determining  the  variation  of  a  given  quantity  between 
neighboring  grid  locations.  With  respect  to  the  variation  of  velocity,  this  problem 
can  be  remedied  by  employing  a  staggered  grid.  The  premise  of  the  staggered  grid 
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is  that  there  may  be  an  advantage  to  storing  some  quantities  at  the  main  grid 
location  while  storing  other  variables  elsewhere.  A  staggered  grid  is  employed  in 
the  FORDC-2  code. 

A  representative  segment  of  the  finite-difference  grid  is  shown  in  Figure  3.1. 
The  grid  is  constructed  in  the  Cartesian  coordinate  system  in  two  dimensions  ( x , 
y ).  The  capital  letters  (E,  W,  N,  and  S)  mark  the  locations  of  the  main  grid 
points  with  respect  to  the  point  of  interest,  point  P.  Note  that  the  symbols  for  the 
cells  neighboring  point  P  are  derived  from  the  points  on  a  compass.  It  is  at  the 
main  grid  points  that  the  properties  of  the  fluid  are  stored,  with  the  exception  of 
the  velocity  components  u  and  v.  The  dashed  lines  represent  control  volume  faces, 
shown  as  lower  case  letters,  which  are  the  boundaries  of  the  control  volume  at  point 
P.  The  u  velocity  components  are  stored  at  the  e  and  w  faces  and  the  v  velocity 
components  are  stored  at  the  n  and  s  faces.  It  is  the  storage  of  the  velocities  at  the 
control  volume  faces,  removed  from  the  main  grid  points,  that  creates  the  staggered 
grid.  The  staggered  grid  can  simply  be  thought  of  as  three  different  grids  being 
employed  on  the  same  calculation  domain:  one  grid  is  employed  for  the  primary 
grid  quantities,  the  second  for  u  velocities,  and  the  third  for  v  velocities. 

The  advantage  to  be  gained  by  this  grid  is  that  the  mass  flow  across  a  control 
volume  face  can  be  calculated  directly.  If  the  velocities  were  stored  at  the  main  grid 
points  with  the  other  variables,  the  values  of  the  velocities  at  the  grid  faces  would 
arise  from  an  interpolation  and  introduce  an  unnecessary  approximation  in  the  cal¬ 
culation  of  the  mass  flow  through  a  control  volume.  Since  satisfaction  of  continuity 
is  fundamental  to  the  solution  of  any  flow  problem,  it  is  critical  to  calculate  accu¬ 
rate  mass  flow  rates.  The  employment  of  the  staggered  grid  aids  in  the  achievement 
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of  that  end.  This  type  of  grid  also  allows  the  direct  calculation  of  the  pressure  gra¬ 
dient  for  a  given  control  volume. 

The  finite-difference  grid  is  constructed  for  the  entire  calculation  domain  by 
the  specification  of  the  control  volume  face  locations.  Once  the  face  locations  have 
been  identified,  the  main  grid  points  are  placed  such  that  they  are  centered  between 
the  control  volume  faces  in  each  direction.  The  spacing  between  the  faces  can  be 
either  uniform  or  nonuniform.  The  advantage  of  the  nonuniform  grid  is  that  the 
grid  can  be  made  more  dense  in  areas  of  interest  without  making  the  entire  grid 
dense.  This  allows  one  to  get  more  detailed  information  about  the  flow  field  in  the 
areas  of  interest  without  increasing  computer  run  time. 


The  uniqueness  of  a  given  problem  is  determined  by  the  specification  of  bound¬ 
ary  conditions.  In  the  case  of  pressure  and  temperature,  the  boundaries  were  typi¬ 
cally  set  to  take  on  ambient  conditions.  However,  the  boundary  conditions  imposed 
on  velocities  are  the  primary  interest  in  this  case. 

The  inlet  and  exit  boundary  conditions  are  critical  to  calculating  accurate  flow 
fields.  In  these  studies,  the  inlet  velocity  plane  is  always  situated  at  the  location 
of  minimum  x,  and  the  outlet  velocity  plane  is  always  situated  at  the  location  of 
maximum  x.  The  inlet  velocity  distribution  is  simply  specified  to  match  the  desired 
velocity  profile  at  the  inlet  of  the  calculation  domain.  For  these  studies,  the  inlet  is 
always  a  uniform  velocity,  U <»•  The  exit  boundary  velocity  must  be  such  that  the 
continuity  equation  is  satisfied.  This  result  can  be  achieved  in  two  manners  by  the 
FORDC-2  code.  First,  the  velocity  at  the  exit  plane  can  be  specified  directly  to 
satisfy  conservation  of  mass.  The  second  scheme  available  adds  some  sophistication 
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to  the  calculation  of  the  exit  velocity  profile,  but  also  carries  a  restriction  on  the 
direction  of  the  exiting  flow. 

The  second  scheme  for  setting  the  exit  boundary  condition  is  intended  to  allow 
the  finite-difference  algorithm  to  calculate  the  proper  exit  velocity  profile.  The 
scheme  involves  calculating  the  mass  flow  across  the  plane  just  ahead  of  the  exit 
plane,  scaling  all  of  the  velocities  on  this  plane  to  match  the  entering  mass  flow,  and 
then  placing  these  scaled  velocities  at  the  exit  plane.  In  this  manner,  the  velocity 
profile  at  the  exit  is  similar  to  the  profile  just  upstream  of  the  exit  and  overall 
continuit}'  is  satisfied.  However,  one  must  be  cautious  to  place  the  exit  boundary 
far  downstream,  for  any  backflow  at  the  exit  will  prevent  this  scheme  from  working 
properly. 

The  other  exterior  boundaries  to  the  flow  are  handled  by  either  external  walls 
or  symmetry  planes.  If  the  boundary  is  considered  to  be  a  wrall,  the  velocity  com¬ 
ponents  are  set  to  zero  at  the  wall.  However,  if  the  boundary  is  placed  at  a  fluid- 
fluid  interface,  then  the  proper  boundary  condition  is  a  symmetry  plane  or  slip  wall. 
The  symmetry  boundary  condition  is  imposed  by  setting  the  finite-difference  coef¬ 
ficient  that  links  the  boundary  point  to  an  interior  cell  equal  to  zero.  This  insures 
that  the  velocity  stored  at  the  symmetry  boundary  has  no  effect  on  the  flow  field. 
This  is  typically  used  in  a  situation  where  the  boundary  could  also  have  been  con¬ 
sidered  to  be  at  freestream  conditions.  For  this  reason,  the  location  of  a  symmetry 
boundary  is  typically  far  from  the  body  about  which  the  fluid  is  flowing.  Of  course, 
if  there  is  an  actual  symmetry  in  the  flow  field,  then  the  symmetry  boundary  can 
be  used  at  a  symmetry  plane  to  shrink  the  calculation  domain.  This  allows  one  to 
achieve  greater  resolution  in  the  flow  field  while  still  maintaining  the  same  number 


of  control  volumes.  A  simple  example  where  an  actual  symmetry  plane  could  be 
used  is  the  flow  between  two  flat  plates,  or  Poiseuille  flow. 

Another  important  consideration,  though  not  technically  a  boundary  condition, 
is  the  internal  wall.  For  these  studies,  internal  walls  are  always  considered  to  be 
impermeable.  The  walls  are  placed  at  control  volume  faces  since  the  specification 
of  the  wall  requires  that  the  velocity  component  perpendicular  to  the  wall  be  set 
to  zero.  The  finite-difference  coefficient  that  links  the  velocities  parallel  to  the 
wall  on  either  side  of  the  wall  is  set  to  zero.  This  prevents  these  velocities  from 
wrongly  influencing  each  other  across  a  wall.  Finally,  the  viscous  effects  that  appear 
near  walls  are  accounted  for  by  employing  a  wall  function  (Nallasamy  1985).  The 
numerical  implementation  of  the  wall  function  is  outlined  by  Carroll  et  al. 

Turbulence  Model 

The  accurate  modelling  of  turbulent  flows  is  a  particularly  challenging  prob¬ 
lem.  Theoretically,  solution  of  the  full  Navier-Stokes  equations  should  provide  an 
accurate  result  with  the  full  effects  of  turbulence  present.  The  problem  lies  in  the 
fact  that  the  scale  of  the  turbulent  fluctuations  are  on  the  order  of  10-3  that  of  the 
primary  flow  (Rodi  1984).  For  a  two-dimensional  flow,  this  would  require  an  in¬ 
crease  in  the  number  of  grid  points  by  an  order  of  106  to  handle  the  details  of  the 
turbulence.  Even  on  the  fastest  of  today’s  supercomputers,  solution  of  such  a  large 
set  of  equations  would  require  prohibitive  amounts  of  memory,  time,  and  money. 

To  account  for  the  effects  of  turbulence  without  attempting  the  solution  of  the 
full  Navier-Stokes  equations  one  employs  empirical  turbulence  models.  The  Navier- 
Stokes  equations  are  converted  into  time  averaged  equations  by  assuming  that  the 
flow  quantities  are  composed  of  a  mean  value  and  a  fluctuating  value 
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<t>  =  4>  +  4>' 


(3-9) 


where 

4>  =  mean  value  of  <f) 

<f>‘  =  fluctuating  component  of  <j>. 

Typically  it  is  the  mean  values  of  the  flow  which  are  of  interest.  This  formulation 
spawns  additional  terms  for  which  extra  relations  must  be  derived.  For  the  momen¬ 
tum  equation,  these  additional  terms  are  known  as  the  Reynolds  stresses.  It  is  the 
task  of  the  turbulence  model  to  express  these  terms  as  functions  of  the  mean  flow 
quantities. 

The  turbulence  model  employed  in  the  FORDC-2  program  is  perhaps  the  most 
popular  turbulence  model  used  in  the  computational  fluid  dynamics  codes,  the  fc- 
e  model.  The  basis  of  this  model  is  that  the  turbulent  flow  can  be  characterized 
by  two  quantities,  the  kinetic  energy  of  turbulence,  k,  and  the  dissipation  rate  of 
turbulence  energy,  e.  The  solution  of  the  differential  equations  for  k  and  e  allows 
the  calculation  of  a  turbulent,  or  eddy  viscosity,  fit.  The  eddy  viscosity  is  added  to 
the  laminar  viscosity  to  yield  an  effective  viscosity, 


Me//  =  M  + 


(3.10) 


It  is  this  variation  in  the  local  viscosity  which  accounts  for  the  presence  of  turbulence 
in  a  flow.  This  model  was  derived  empirically,  and  expectedly  does  not  work  well 
for  all  flow  situations.  For  a  more  detailed  description  of  the  fc-e  turbulence  model 
refer  to  Rodi  or  Nallasamy. 
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CHAPTER  IV 


DETERMINATION  OF  AERODYNAMIC 
FORCE  COEFFICIENTS  AND 
SHEDDING  FREQUENCY 

The  following  discussion  outlines  the  methods  for  calculating  various  parame¬ 
ters  associated  with  a  given  flow.  The  methodology  for  calculating  the  drag  coeffi¬ 
cient  C £> ,  the  lift  coefficient  Cl,  and  the  Strouhal  number  St  are  discussed. 

Drag  Coefficient 

In  the  most  general  treatment,  the  aerodynamic  forces  experienced  by  a  body 
are  the  result  of  the  distribution  of  shear  stress  and  pressure  on  the  surface  of  the 
body.  The  drag  force  is  the  force  experienced  by  the  body  in  a  direction  parallel 
to  the  freestream  velocity.  For  the  case  of  a  flat  plate,  if  one  considers  the  plate  to 
be  infinitesimally  thin  the  effects  of  shear  and  pressure  on  the  width  of  the  plate 
can  be  neglected.  Therefore,  the  drag  on  a  plate  which  is  normal  to  the  freestream 
direction  is  simply  a  function  of  the  pressure  distribution  on  the  plate. 

Referring  to  Figure  4.1  one  can  immediately  write  the  drag  force  as 

/  i 

D  =  J pf(y)dA  -  J  pb(y)dA  .  (4.1) 


This  calculation  of  drag  is  dependent  on  a  knowledge  of  the  continuous  distribution 
of  the  pressure  force  over  the  surface  of  the  plate.  However,  in  the  computational 
scheme  the  surface  of  the  plate  has  been  broken  into  a  discrete  number  of  elements, 
N.  An  example  of  the  discretized  pressure  distribution  is  given  in  Figure  4.2. 
Defining  the  depth  of  the  plate  as  a  unit  1,  the  differential  area  then  becomes 

dA  =  dy  ■  1  =  dy  .  (4-2) 

One  can  then  write  the  drag  in  terms  of  the  discretized  pressures  as 

N 

D  ~  J^(P/  -  Pb)i&Vi  •  (4.3) 

i=l 

This  expression  can  be  simplified  by  denoting  the  pressure  differential  as 

&Pi  =  {Pf~Pb)i  (4.4) 

which  gives  an  expression  for  the  drag 

N 

D=^ApiAyi  .  (4.5) 

i=i 

This  expression  can  be  further  simplified  by  noting  that  the  plate  is  discretized  into 
elements  of  equal  length  Ayl.  This  yields  the  simplified  drag  expression 

N 

D  =  Ay  ^  A  pi  .  (4.6) 

i=i 

The  task  is  then  to  cast  the  drag  force  in  terms  of  a  drag  coefficient.  Recall 
that  the  definition  of  the  drag  coefficient  is 


Substituting  the  expressions  for  the  drag  and  the  area  yields 


oriented  normal  to  the  freestream  direction,  the  lift  is  due  only  to  the  shear  forces 
which  act  on  the  surface  of  the  fiat  plate. 

Figure  4.3  shows  the  distribution  of  the  shear  stress  r{y)  on  the  surface  of  a 
flat  plate.  The  shear  stress  for  a  Newtonian  fluid  in  the  coordinate  system  defined 
in  Figure  4.3  is 


(4.13) 


The  velocity  gradient  dv / dx  can  be  computed  by  simply  dividing  the  v  velocity  that 
neighbors  the  surface  by  the  distance  between  the  velocity  location  and  the  surface. 
The  gradient  may  be  calculated  in  this  manner  because  the  velocity  at  the  surface 
must  be  zero. 

Turning  now  to  the  calculation  of  the  lift  force,  one  can  immediately  write  that 


l  l 

L  =  J  Ty,(y)dA+  J  Tyb{y)dA  ■  (4-14) 

0  0 

Proceeding  in  the  same  manner  as  that  used  in  the  calculation  of  Cp,  the  lift  for 
equal  Ayt  can  be  written  as 


N 

L  =  Ay^(rv/  ^ryJ,  .  (4.15) 

t=i 

Recalling  the  definition  of  the  lift  coefficient  as 


(4-16) 


the  lift  and  area  expressions  can  be  substituted  to  yield 


i  , 


Aerodynamic  Forces  on  a  Rectangular  Cylinder 


With  the  lift  and  drag  on  a  flat  plate  resolved,  the  task  now  is  to  extend 


those  forces  to  a  rectangular  cylinder.  The  pressure  and  shear  forces  which  act  on 


the  cylinder  are  shown  in  Figure  4.5.  Note  that  for  the  rectangular  cylinder  lift 


and  drag  both  have  force  components  due  to  pressure  and  shear  stress.  Recalling 


the  definitions  of  Cp  and  Cl  once  more,  and  using  the  forms  of  the  lift  and  drag 


forces  given  in  the  previous  section,  the  lift  and  drag  coefficients  on  the  rectangular 


cylinder  are  easily  derived. 


Drag  on  the  cylinder  oriented  normal  to  the  freestream  velocity  is  due  to  the 


pressure  on  the  front  and  back  surfaces  and  the  shear  stresses  on  the  upper  and 


lower  surfaces.  Once  again  the  spacing  of  the  Ax{  and  Ay;  terms  is  uniform,  with 


N  elements  composing  a  vertical  surface  and  M  elements  composing  a  horizontal 


surface.  An  expression  for  the  drag  in  terms  of  average  pressures  and  stresses  can 


be  written  immediately  as 


D  —  N(pf  -  pb)Ay  -f  M(tXu  +  tXi  )Ai  . 


(4.23) 


Recalling  that  N Ay  is  l  and  M Ax  is  c,  the  drag  coefficient  becomes 


^  (Pf  ~  Pb )  ,  (TIu  +  TXI)C 
Od  = - 1 - ; - 

9oo  9oo‘ 


(4.24) 


Lift  on  the  cylinder  is  due  to  the  pressure  on  the  upper  and  lower  surfaces  and 


the  shear  stresses  on  the  front  and  back  surfaces.  The  lift  can  then  be  written  as 


L  =  M{pi  -  pu) Ax  +  N(Tyf  +  ryjAy  . 


'sV/k 
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The  lift  coefficient  is  then 

~  Pu)c  ^  (ry f  Tyh 

Qool  Qoo 

The  expressions  for  calculating  Cp  and  Cl  for  both  a  flat  plate  and  a  rectangular 
cylinder  perpendicular  to  the  freestream  flow  have  been  determined.  These  equa¬ 
tions  are  the  relations  programmed  into  the  FORDC-2  code  for  the  determination 
of  the  aerodynamic  force  coefficients. 

Shedding  Frequency 

The  method  for  determining  the  shedding  frequency  is  derived  from  the  work  by 
Younis  (1988).  He  proposed  to  determine  the  shedding  frequency  by  examining  the 
periodic  fluctuations  in  three  fundamental  quantities.  He  determined  the  frequency 
n  by  examining 

1.  the  axial  velocity  at  two  locations  behind  the  body  as  shown  in  Figure  4.6,  the 
velocities  being  denoted  by  u\  and  1x2, 

2.  the  drag  coefficient  C p,  and 

3.  the  lift  coefficient  Cl- 

For  the  slightly  more  complex  case  of  a  rectangular  cylinder,  the  axial  velocity 
is  measured  at  four  locations  as  shown  in  Figure  4.7  because  of  the  possibility  of 
vortices  being  shed  from  either  the  front  or  rear  of  the  cylinder.  Note  that  the  flow 
field  must  be  allowed  to  reach  a  cyclic  state  before  the  frequency  can  be  determined 
by  this  method. 

The  shedding  frequency  can  be  extracted  directly  from  plots  of  either  u  1,  u2, 
Cp ,  or  Cl  versus  time.  This  method,  with  one  exception,  assumes  that  vortices 
are  shed  at  the  same  frequency  at  which  the  aforementioned  quantities  fluctuate. 


(4.26) 


Knowledge  of  the  shedding  frequency  immediately  leads  to  the  Strouhal  number  by 
the  relation 

St  =  ~  .  (4.27) 

V  oo 

Younis  notes  that  the  drag  oscillates  at  twice  the  frequency  of  the  lift,  which 
oscillates  at  the  shedding  frequency.  Plots  of  the  fluctuations  in  Cd  and  C l  provided 
by  Younis  are  given  as  Figure  4.8  and  clearly  demonstrate  the  higher  frequency 
oscillation  in  the  drag.  Similar  plots  obtained  in  this  research  for  a  flat  plate  at 
Re  —  7826  are  given  in  Figures  4.9-4.12,  and  the  higher  frequency  oscillation  of  the 
drag  is  again  evident. 


FIGURE  4.5:  Dis  !  ation  of  Pressure  and  Shear  Stress 
on  a  Rectangular  Cylinder. 


FIGURE  4.6:  Location?  for  Transient  Velocity 
Measurements  for  a  Flat  Plate. 


FIGURE  4.11:  Variation  in  Drag  Coefficient  Over  Time 


CHAPTER  V 


RESULTS  FOR  VORTEX  SHEDDING 
FROM  A  FLAT  PLATE 


The  subsequent  discussion  describes  the  results  obtained  for  the  case  of  vortex 
shedding  from  a  flat  plate  in  a  uniform  flow.  Several  comparisons  are  made  with 
data  obtained  by  Roshko  for  this  geometrical  configuration.  The  effects  of  changing 
parameters  associated  with  the  finite-difference  solution  of  the  flow  field  are  also 
examined. 


Computational  Setup 

The  numerical  prediction  of  vortex  shedding  from  a  flat  plate  in  a  uniform  flow 
was  achieved  on  for  different  grid  configurations.  The  coarser  of  the  grids  used  in 
the  study  was  a  60  x  50  uniform  grid  as  shown  in  Figure  5.1.  Overall  physical 
dimensions  of  the  calculation  domain  were  0.06  m  x  0.03  m.  The  flat  plate  was 
composed  of  10  vertical  cells,  had  a  length  of  0.006  m,  and  was  centered  vertically 
within  the  calculation  domain.  The  spacing  between  the  plate  and  the  upper  and 
lower  boundaries  was  approximately  two  plate  lengths.  An  inlet  boundary  was 
provided  at  approximately  2.5  plate  lengths  upstream  of  the  plate,  and  the  exit 
boundary  was  placed  approximately  7.5  plate  lengths  downstream  of  the  plate. 
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A  finer  grid  was  also  employed  in  the  solution  of  this  problem.  The  second  grid 
was  a  highly  nonuniform  150  x  80  grid  as  shown  in  Figure  5.2.  Nonuniformity  was 
introduced  to  concentrate  the  control  volumes  in  the  area  of  greatest  interest,  that 
area  being  the  region  directly  downstream  of  the  plate.  Also,  the  nonuniformity 
allows  for  the  placement  of  relatively  large  cells  near  the  boundaries  which  locate 
the  boundaries  farther  away  from  the  plate  to  reduce  the  influence  of  boundary 
conditions.  The  overall  dimensions  of  the  calculation  domain  were  0.175  m  x  0.140 
m.  The  flat  plate  was  composed  of  16  equally  spaced  vertical  cells  which  gave  a  plate 
length  of  0.0102  m.  Vertical  spacing  between  the  plate  and  the  outer  boundaries 
was  approximately  6.4  plate  lengths  at  both  the  top  and  bottom.  Upstream  of  the 
plate  by  about  two  plate  lengths  was  an  inlet  boundary,  and  approximately  15  plate 
lengths  downstream  of  the  plate  was  an  exit  boundary. 

Boundary  conditions  for  the  two  grid  sizes  were  similar  in  nature.  The  upper 
and  lower  boundaries  were  assigned  a  symmetry  (slip  wall)  condition.  The  inlet 
boundary  was  a  uniform  flow  directed  normal  to  the  plate  with  the  magnitude  of  the 
velocity  being  varied  from  0.3  m/s  to  3.0  m/s  to  achieve  a  range  of  Reynolds  numbers 
of  approximately  2000  <  Re  <  30000.  The  exit  boundary  was  set  by  allowing 
the  boundary  to  assume  a  velocity  profile  similar  to  the  one  directly  upstream  as 
described  in  Chapter  III. 

The  FORDC-2  program  was  executed  on  the  VAX  8650  and  VAX  11/780 
computers  available  at  the  Texas  Tech  computing  facilities.  The  problem  was  solved 
isothermally  with  the  calculated  variables  being  u,  v,  and  p.  Runs  were  also  made 
for  the  coarse  grid  case  with  the  k-e  turbulence  model  included  in  the  solution, 
which  results  in  the  addition  of  the  variables  k  and  e  to  the  process.  Turbulence 


was  employed  to  determine  whether  there  was  a  significant  difference  in  the  results 
between  the  turbulent  and  laminar  cases.  The  different  grids  were  employed  to 
examine  the  effects  of  a  refined  grid  on  the  solution  and  to  determine  whether  the 
extra  effort  required  to  create  and  run  the  more  complex  grid  was  justified.  Finally, 
the  time  step,  At,  was  adjusted  on  a  case  by  case  basis  to  ensure  convergence  of 
the  transient  solution. 

Method  of  Inducing  Vortex  Shedding 

A  significant  problem  arose  in  computationally  simulating  the  Karman  vortex 
street.  Geometrically,  the  problem  has  a  symmetry  about  the  horizontal  line  which 
bisects  the  flat  plate.  Since  the  finite-difference  model  is  simply  a  mathematical 
model,  the  flow  field  that  results  would  share  this  symmetry  about  the  horizontal 
bisector.  This  symmetry  would  seemingly  exclude  the  vortex  street  pattern  from 
the  set  of  possible  solutions. 

Fortunately,  this  problem  can  be  overcome  by  removing  the  symmetry  in  the 
initial  conditions.  To  make  the  flow  initially  asymmetric,  a  perturbation  was  added 
to  the  flow  in  the  form  of  a  momentum  source.  The  momentum  source  was  included 
for  the  first  three  time  steps  and  was  located  at  a  control  volume  directly  behind 
and  at  the  upper  limit  of  the  plate.  The  location  of  the  source  was  approximately 
the  same  location  at  which  the  transient  velocity,  u\ ,  was  measured  as  shown  in 
Figure  4.6.  Inclusion  of  this  momentum  source  proved  to  be  an  adequate  method 
of  introducing  asymmetry  to  the  flow;  consequently,  the  vortex  street  was  predicted 
by  the  program. 
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Qualitative  Analysis  of  Transient 
Vector  Fields 

Prediction  of  the  Karman  vortex  street  by  the  FORDC-2  program  is  clearly 
demonstrated  by  Figures  5. 3-5. 7.  These  figures  show  the  absolute  velocity  vectors 
for  equally  separated  increments  in  time.  The  five  plots  progress  by  steps  of  a 
quarter  of  the  vortex  shedding  period,  P,  to  show  the  progression  of  the  vector 
field  with  time.  Absolute  velocity  vectors  are  displayed  simply  because  the  vortex 
shedding  pattern  is  far  more  apparent  when  shown  in  this  manner.  The  velocity 
vectors  are  transformed  by  the  simple  relation 

V»bs  —  V  Uqo  .  (5.1) 

At  t  =  0  (Figure  5.3),  the  familiar  alternating  vortex  pattern  of  the  vortex 
street  is  clearly  visible.  In  fact,  the  pattern  shown  has  much  the  same  character 
as  that  given  in  Figure  1.1.  There  is  a  pronounced  spreading  of  the  vortices  in  the 
vertical  direction  as  they  move  downstream.  Also  note  that  the  vortices  seemingly 
grow  in  size  as  they  progress  downstream.  These  characteristics  indicate  that  the 
finite-difference  code  accurately  predicts  the  vortex  shedding  pattern. 

Progressing  to  t  =  1/4P  (Figure  5.4\  notice  that  the  entire  flow  pattern  from 
t  =  0  has  moved  downstream.  There  is  also  a  disturbance  at  the  top  of  the  plate 
which  is  the  beginning  of  the  formation  of  the  next  vortex.  At  t  =  1/2 P  (Figure 
5.5),  the  vortex  on  the  top  of  the  plate  has  formed  and  is  clearly  recognizable.  As 
expected,  ct i.  I  =  1/2 P  the  vector  field  appears  to  be  a  mirror  image  about  the 


horizontal  bisector  of  the  field  at  t  =  0.  The  vector  field  at  t  =  3/4 P  (Figure  5.6) 
indicates  the  beginning  of  the  formation  of  a  new  vortex  at  the  bottom  of  the  plate. 


wvtvrcw  ww  v^vicxx; 


|L!< 


Finall}',  at  t  =  P  (Figure  5.7)  the  vector  field  appears  to  be  the  same  field  calculated 
at  t  =  0.  Therefore,  the  vector  field  has  cyclicly  progressed  back  to  its  initial  state, 
and  the  predicted  vortex  shedding  is  indeed  periodic. 

Comparison  of  Numerical  and  Experimental 


Strouhal  Numbers 

Among  the  most  critical  parameters  associated  with  periodic  vortex  shedding 
is  the  Strouhal  number.  Roshko  provided  extensive  data  on  the  Strouhal  numbers 
for  different  geometries,  including  the  flat  plate.  Figure  5.8  offers  a  quantitative 
comparison  of  the  experimental  results  of  Roshko  with  the  results  from  this  study. 
The  data  presented  graphically  are  also  given  in  tabular  form  as  Tables  5.1  and  5.2. 

A  quick  examination  of  Figure  5.8  indicates  that  the  numerical  model  predicted 
Strouhal  numbers  comparable  to  those  measured  by  Roshko.  However,  a  closer 
examination  reveals  some  important  trends  in  the  data.  For  the  detailed  150  x  80 
grid,  the  prediction  of  Strouhal  number  was  quite  accurate  up  to  a  Reynolds  number 
of  approximately  104.  For  Reynolds  numbers  greater  than  104,  the  model  seemed 
to  lose  accuracy  as  the  predicted  Strouhal  number  decreased.  In  any  case,  the  error 
between  the  predicted  and  experimental  results  never  exceeded  15  percent. 

For  the  case  of  the  60  x  50  grid  the  same  trend  was  not  apparent.  For  the 
entire  Reynolds  number  range,  the  coarser  grid  appears  to  have  predicted  Strouhal 
numbers  more  consistent  with  Roshko’s  results.  The  Strouhal  numbers  predicted 
on  the  coarse  grid  were  generally  greater  than  those  predicted  on  the  fine  grid. 
Addition  of  the  k-e  turbulence  model  generally  caused  a  slight  increase  in  the 
Strouhal  number;  however,  that  trend  was  not  substantiated  for  all  cases. 


Overall,  the  results  were  favorable,  especial!}'  when  compared  with  the  work 
done  by  Younis.  In  a  similar  study,  Younis  predicted  Strouhal  numbers  for  a  flat 
plate  in  a  range  of  0.165-0.185.  Younis’  values  would  appear  to  be  an  overestimate 
of  the  Strouhal  number  when  compared  with  the  experimental  data  provided  by 
Roshko.  Overestimation  of  the  Strouhal  was  not  a  problem  encountered  in  the 
current  study. 


Wake  Pressure  Profiles 


Roshko  provided  experimental  results  for  the  time  averaged  pressure  distribu¬ 
tion  on  the  centerline  of  the  wake.  His  data  was  for  the  specific  Reynolds  number 
of  Re  =  14500,  and  the  data  was  given  in  terms  of  the  nondimensional  quantities 
Cp  and  x/l.  The  quantity  x  is  measured  from  the  plate  with  positive  x  extending 
downstream.  A  comparison  with  Roshko’s  data  was  made  for  both  grid  sizes  and 
for  the  case  where  turbulence  was  included. 

Figure  5.9  gives  the  variation  in  Cp  with  downstream  distance  for  the  fine  grid. 
The  most  obvious  difference  between  the  computational  and  experimental  curves  is 
the  magnitude  of  the  suction  predicted  by  the  numerical  model  just  downstream  of 
the  plate.  The  numerical  model  predicted  a  base  pressure  coefficient  nearly  twice  as 
large  as  Roshko  measured.  Note  that  the  numerical  model  predicted  the  minimum 
pressure  to  occur  at  x/l  ~  0.7  where  the  experimental  results  indicate  that  the 
minimum  pressure  occurs  at  r/1  ~  2.0.  Minimum  pressure  values  were  also  quite 
different,  with  the  predicted  minimum  taking  a  value  approximately  twice  that  of 
the  experimental  result.  However,  this  does  not  necessarily  indicate  an  error  in  the 
numerical  model,  because  Roshko  cited  difficulties  in  obtaining  the  wake  pressure 


measurements.  Furthermore,  note  that  the  agreement  between  the  curves  for  x/l  > 
2.0  was  quite  good. 

The  anomalies  cited  for  the  fine  grid  solution  are  also  apparent  in  the  coarse 
grid  solution  as  shown  in  Figure  5.10.  The  predicted  base  pressure  coefficients 
for  the  two  grids  are  almost  identical.  Minimum  pressure  for  the  two  cases  also 
occurred  at  nearly  the  same  value  of  x/l,  though  the  coarse  grid  solution  did  not 
predict  as  low  a  value  of  Cp  as  did  the  fine  grid  solution.  The  primary  difference 
between  the  solutions  for  the  two  grid  lies  in  the  accuracy  of  the  solution  for  values  of 
x/l  >  2.0.  For  this  region,  the  fine  grid  solution  was  much  closer  to  the  experimental 
solution  than  was  the  coarse  grid  solution.  Nonetheless,  the  coarse  grid  solution  did 
approximate  this  part  of  the  curve  with  reasonable  accuracy. 

The  addition  of  turbulence  to  the  coarse  grid  solution  appears  to  make  little 
difference  in  the  solution.  The  pressure  profile  for  the  turbulent  case  (Figure  5.11) 
is  virtually  identical  to  the  laminar  case.  This  result  indicates  that  the  addition  of 
turbulence  has  little  effect  on  the  prediction  of  wake  pressures. 

Comparison  of  Numerical  and  Experimental 
Aerodynamic  Force  Coefficients 

Measurement  of  the  time  averaged  lift  and  drag  coefficients  was  carried  out 
for  each  case  run.  As  expected,  the  average  lift  coefficient,  Ciavg,  was  zero  for  all 
runs  though  instantaneous  values  of  lift  were  nonzero.  The  zero  value  of  lift  was 
expected  since  this  case  is  for  a  flat  plate  at  zero  angle  of  attack. 

A  graphical  display  of  the  calculated  average  drag  coefficients  is  given  as  Figure 
5.12.  Also  included  on  the  plot  are  drag  coefficients  derived  from  Roshko’s  data 
and  the  equations  given  in  Chapter  II.  Note  the  extreme  difference  between  the 
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numerical  solution  and  the  experimental/theoretical  values  provided  by  Roshko. 
The  numerical  values  of  Co  exceed  Roshko’s  values  by  more  than  50  percent,  and 
this  may  well  be  tied  to  the  much  lower  base  pressure  predicted  by  the  computations. 
However,  also  note  that  the  predicted  values  agree  favorably  with  the  numerical 
results  of  Younis  who  found  Co' s  in  the  range  of  3.25-3.65. 

From  a  numerical  standpoint  there  were  also  significant  differences  in  the 
predicted  drag  coefficients.  The  fine  grid  solution  estimated  drag  coefficients  on  the 
order  of  10  percent  greater  than  the  coarse  grid.  This  fact  may  be  related  to  the 
lower  minimum  pressure  found  in  the  wake  for  the  fine  grid  solution.  Turbulence 
appeared  to  have  little  effect  on  the  prediction  of  drag,  though  its  addition  did  seem 
to  generally  increase  the  drag  coefficient  by  about  1  percent. 

Effect  of  Iterative  Time  Step  Size 


From  a  numerical  standpoint,  the  selection  of  the  iterative  time  step  At  can 
be  critical.  Figure  5.13  presents  the  variation  in  the  predicted  Strouhal  number  for 
various  iterative  time  step  sizes.  The  horizontal  line  indicates  the  experimentally 
measured  Strouhal  number  for  this  particular  Reynolds  number,  Re  =  14500.  A 
best  fit  linear  approximation  of  the  data  is  also  displayed.  Note  that  the  curvature 
in  the  curve  fit  is  due  to  the  logarithmic  axis. 

From  this  figure  one  can  clearly  see  the  dependence  of  the  solution  on  the  time 
step.  The  prediction  of  the  Strouhal  number  deteriorates  markedly  as  the  time  step 
is  increased.  For  this  particular  case,  any  further  increase  in  the  time  step  above  A l 
=  .0008  led  to  divergence  of  the  solution.  Also  note  that  as  A t  — *  0  the  predicted 
solution  asymptotically  approaches  the  experimental  solution.  These  indicate  that 
the  choice  of  a  sufficiently  small  value  of  At  is  critical  in  obtaining  accurate  solutions. 
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FIGURE  5.8:  Strouhal  Number  vs.  Reynolds  Number 
for  a  Flat  Plate. 


TABLE  5.1:  Numerical  Shedding  Frequency  and 
Aerodynamic  Force  Data  for  the 
150  x  80  Grid. 
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FIGURE  5.9:  Centerline  Wake  Pressure  Distribution  for 


FIGURE  5.10:  Centerline  Wake  Pressure  Distribution  for 


FIGURE  5.11:  Centerline 


FIGURE  5.13:  Effect  of  Iterative  Time  Ste] 
Predicted  Strouhal  Number 


CHAPTER  VI 


RESULTS  FOR  VORTEX  SHEDDING  FROM 
A  FLAT  PLATE  NEAR  A  GROUND  PLANE 

The  following  discussion  describes  the  results  obtained  for  the  case  of  a  flat 
plate  near  a  ground  plane.  A  comparison  is  made  between  the  predicted  Strouhal 
numbers  and  the  function  obtained  by  Matty.  The  effect  of  the  presence  of  the 
ground  plane  on  the  aerodynamic  force  coefficients  is  also  examined. 

Computational  Setup 

The  numerical  prediction  of  vortex  shedding  from  a  flat  plate  near  a  ground 
plane  was  carried  out  on  only  one  grid  configuration.  This  is  due  to  the  fact  that 
results  in  Chapter  V  indicated  that  refinement  of  the  grid  yielded  only  limited 
improvements  in  the  solution.  However,  the  60  x  50  grid  used  in  the  previous  case 
proved  to  be  unacceptable  due  to  a  lack  of  spacing  between  the  plate  and  the  exit 
boundary.  To  remedy  this  problem,  cells  were  added  to  the  grid  downstream  of  the 
plate. 

The  grid  employed  for  this  problem  was  a  80  x  50  grid  as  shown  in  Figure  6.1 
with  overall  dimensions  of  0.08  m  x  0.03  m.  The  plate  was  composed  of  12  vertical 
cells  and  had  a  length  of  0.0072  m.  The  inlet  boundary  was  located  approximately 
two  plate  heights  upstream  of  the  plate,  and  the  exit  was  located  approximately 


nine  plate  heights  downstream  of  the  plate.  The  location  of  the  bottom  boundary 
relative  to  the  plate  was  governed  by  the  parameter  A  which  was  assigned  the  values 
of  0,  1/4,  1/2,  and  1.  The  spacing  between- the  plate  and  the  upper  boundary  varied 
between  2  and  3  plate  heights  for  different  values  of  A. 

Boundary  conditions  for  this  case  were  similar  to  those  presented  in  Chapter 
V.  The  primary  difference  was  the  placement  of  an  impermeable  wall  at  the  bottom 
boundary.  Otherwise,  the  boundaries  remained  the  same  w’ith  the  exception  of  the 
magnitude  of  the  inlet  velocities.  The  inlet  velocity  was  varied  between  2.0  m/s  and 
20.0  m/s  to  give  a  range  of  Reynolds  numbers  of  approximately  1000  <  Re  <  15000. 
This  range  of  Reynolds  numbers  was  chosen  to  match  the  range  examined  by  Matty. 

Vortex  shedding  was  induced  in  the  same  method  explained  in  Chapter  V. 
Turbulence  was  not  included  in  this  case  due  to  the  previous  results  which  indicated 
that  the  addition  of  turbulence  was  of  little  consequence.  Again,  the  solution  was 
isothermal,  and  the  calculated  variables  were  u,  v,  and  p. 


ualitative  Analysis  of  the  Flow  Field 


Inducing  vortex  shedding  in  the  wake  of  a  flat  plate  near  a  ground  plane  was 
not  successful  for  all  cases.  For  the  case  of  A  =  0,  vortex  shedding  was  never 
achieved.  To  achieve  vortex  shedding,  the  plate  had  to  be  moved  to  a  distance  of 
one  plate  length  from  the  wall.  At  this  distance  from  the  wall,  vortex  shedding 
occurred  readily.  To  achieve  shedding  for  other  values  of  A,  the  velocity  field  of  the 
A  =  1  case  was  imposed  as  an  initial  condition. 

Figure  6.2  and  A. 2  respectively  present  representative  plots  of  the  velocity 
vectors  and  streamlines  for  an  instant  in  time.  This  particular  plot  is  for  the  case 
of  A  =  1.  The  main  difference  to  note  between  this  case  and  the  cases  presented  in 
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Chapter  V  is  a  decrease  in  the  wake  width.  This  is  apparently  due  to  the  influence  of 
the  wall  which  deters  the  vortices  shed  from  the  bottom  of  the  plate  from  spreading. 

Comparison  of  Numerical  and  Experimental 


Strouhal  Numbers 

Data  obtained  for  the  various  values  of  A  are  provided  as  Tables  6. 1-6. 3.  A 
comparison  of  the  predicted  Strouhal  number  data  versus  the  function  St(Re.X) 
obtained  by  Matty  is  provided  in  Figures  6. 3-6. 5. 

For  the  greatest  separation  between  the  plate  and  the  ground  plane  at  A  =  1 
(Figure  6.3),  the  numerical  prediction  of  the  Strouhal  number  for  the  entire  range 
of  Reynolds  numbers  was  fairly  accurate.  As  the  Reynolds  number  was  increased, 
the  numerical  solution  tended  to  give  a  slightly  higher  value  of  the  Strouhal  number 
than  the  function  St(Re,  A).  However,  the  Strouhal  number  generally  decreased  for 
increasing  Reynolds  number  as  indicated  by  Matty’s  study. 

As  the  value  of  A  was  decreased  to  A  =  1/2  (Figure  6.4),  the  trend  changed 
considerably.  Predicted  values  of  the  Strouhal  number  were  lower  than  those  given 
by  the  function  St(Re,  A)  for  low  values  of  the  Reynolds  number.  However,  as  the 
Reynolds  number  was  increased  the  numerical  solution  was  quite  accurate.  The 
trend  of  decreasing  Strouhal  number  with  increasing  Reynolds  number  was  also 
supported  by  the  prediction.  The  results  were  also  similar  for  A  =  1/4  (Figure  6.5), 
though  predicted  Strouhal  numbers  slightly  underestimated  the  empirical  relation. 
Nonetheless,  the  trend  of  decreasing  Strouhal  number  with  increasing  Reynolds 
number  was  again  supported. 

The  variation  in  the  Strouhal  number  with  A  for  a  fixed  value  of  the  Reynolds 
number  was  also  examined.  Strouhal  numbers  for  A  =  1/4  were  generally  greater 
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than  those  for  A  =  1/2.  Matty’s  data  indicates  a  similar  behavior  between  these  val¬ 
ues  of  A.  However,  the  predicted  Strouhal  numbers  at  A  =  1  were  significantly  higher 
than  those  at  A  =  1/2,  though  Matty’s  study  indicates  that  these  values  should  be 
nearly  equal.  The  difference  here  may  well  lie  in  the  fact  that  Matty’s  study  was 
for  a  three-dimensional  plate,  while  the  numerical  study  was  two-dimensional. 

Variation  in  Aerodynamic  Forces 
with  Ground  Clearance 

A  variation  in  the  aerodynamic  forces  was  noted  for  changes  in  the  parameter 
A.  However,  for  a  given  value  of  A,  there  was  little  apparent  dependence  of  the 
drag  force  on  the  Reynolds  number.  Figure  6.6  presents  the  variation  in  the  drag 
coefficient  with  A.  The  trend  for  the  drag  to  increase  as  A  decreases  was  clearly 
apparent.  At  A  =  0,  the  drag  had  attained  a  value  approximately  50  percent  greater 
than  that  obtained  for  a  flat  plate  removed  from  a  ground  plane.  This  behavior  was 
expected,  for  the  increase  in  the  drag  on  a  body  as  it  approaches  a  ground  plane  is 
well  documented. 

An  increase  in  lift  also  occurred  for  decreasing  ground  heights,  though  the  trend 
was  not  as  prominent  as  it  was  for  the  drag.  Figure  6.7  presents  the  variation  in  the 
lift  coefficient  with  A  revealing  this  relationship.  The  lift  appeared  to  have  a  slight 
dependence  on  the  Reynolds  number,  with  an  increase  in  Reynolds  number  being 
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TABLE  6.2:  JNumerical  Shedding  Frequency  and 

Aerodynamic  Force  Data  for  )  =  1/2. 
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TABLE  6.3:  Numerical  Shedding  Frequency  and 

Aerodynamic  Force  Data  for  A  =  1/4. 
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FIGURE  6.6:  Variation  in  the  Drag  Coefficient 
Ground  Clearance. 


CHAPTER  VII 


RESULTS  FOR  VORTEX  SHEDDING  FROM 
A  RECTANGULAR  CYLINDER 

The  following  discussion  reviews  results  of  the  computational  prediction  of  vor¬ 
tex  shedding  from  rectangular  cylinders.  Computed  Strouhal  numbers  are  com¬ 
pared  to  experimental  values  provided  by  Okajima.  Trends  in  the  variations  of 
aerodynamic  forces  are  also  explored. 

Computational  Setup 

The  numerical  prediction  of  vortex  shedding  from  rectangular  cylinders  was 
achieved  on  the  same  basic  80  x  50  grid  used  for  the  study  of  a  flat  plate  near  a 
ground  plane.  The  primary  differences  between  the  two  grids  were  the  placement 
of  the  body  in  the  calculation  domain  and  the  specification  of  boundary  conditions. 

Figure  7.1  presents  the  finite-difference  grid  for  the  specific  case  of  a  rectangular 
cylinder  with  0  =  3.  The  study  was  executed  for  cylinders  of  varying  length  to  height 
ratios,  with  the  0  taking  on  the  values  1,  2,  and  3.  The  height  of  the  cylinder  was 
composed  of  10  cells  giving  a  height  of  0.006  m.  The  length  of  the  box  was  composed 
of  either  6,  12.  or  18  cells  depending  on  the  respective  value  of  0 ,  either  1,  2,  or 
3.  The  inlet  boundary  was  located  approximately  two  cylinder  heights  upstream  of 
the  cylinder,  and  the  exit  varied  from  4.5-  6.5  cylinders  heights  downstream  of  the 


cylinder  depending  on  the  value  of  j3.  Spacing  between  the  cylinder  and  the  top 
and  bottom  boundaries  was  fixed  at  approximately  two  cylinder  heights. 


Boundary  conditions  for  this  case  were  identical  to  those  for  a  fiat  plate  in  a 
uniform  stream  (Chapter  V),  with  the  exception  of  modified  inlet  velocities.  Inlet 
velocities  were  varied  from  approximately  0.04  m/s  to  2.80  m/s  yielding  a  Reynolds 
number  range  of  200  <  Re  <  15000.  This  range  roughly  corresponds  to  the  range 
of  the  data  presented  by  Okajima. 

Again,  the  solution  was  isothermal,  and  the  k-e  turbulence  model  was  excluded 
from  the  solution  process.  This  provided  a  solution  in  terms  of  the  variables  u,  v. 
and  p. 

Qualitative  Analysis  of  the  Flow  Field 

The  unusual  behavior  exhibited  in  plots  of  Strouhal  number  versus  Reynolds 
number  for  cylinders  with  aspect  rations  in  the  range  of  2  <  /3  <  3  was  theorized 
to  be  dependent  on  an  overall  change  in  the  flow"  field.  It  was  hoped  that  the  nu¬ 
merical  solution  of  this  problem  would  demonstrate  this  change  in  the  flow  geome¬ 
try.  However,  the  change  in  the  flow  pattern  was  not  at  all  evident  in  plots  of  the 
velocity  vectors  and  streamlines  for  the  tested  cylinders. 

Plots  of  the  flow  fields  for  each  of  the  values  of  /3  are  provided  as  Figures  7. 2-7. 4 
and  Figures  A.3-A.5.  The  plots  are  all  quite  similar,  and  show  a  pattern  very  much 
like  that  obtained  behind  a  flat  plate.  These  figures  are  all  representative  of  the 
flow  fields  obtained  over  the  entire  range  of  Reynolds  numbers  for  their  respective 
values  of  (3.  Unfortunately,  the  resolution  of  the  solution  was  not  such  that  the 
region  of  reattachment  along  the  cylinder  length  could  be  seen  for  the  two  latter 
cases  where  f. 3  took  on  the  values  of  2  and  3. 


Strouhal  Numbers 


The  most  interesting  aspect  of  the  shedding  of  vortices  from  a  rectangular 
cylinder  was  the  discontinuity  that  occurred  in  plots  of  the  Strouhal  number  versus 
the  Reynolds  number.  This  behavior  was  predicted  by  the  FORDC-2  code  for  the 
expected  values  of  0.  Data  for  these  cases  is  provided  in  Tables  7. 1-7.3. 

Figure  7.5  presents  the  variation  in  the  Strouhal  number  with  Reynolds  num¬ 
ber  for  a  cylinder  with  0  =  1.  The  plot  indicates  that  the  Strouhal  number  was  rel¬ 
atively  insensitive  to  the  Reynolds  number  over  the  entire  range  of  Reynolds  num¬ 
bers.  Maximum  values  of  the  Strouhal  number  occurred  at  minimum  values  of  the 
Reynolds  number.  This  prediction  is  supported  by  the  data  provided  by  Okajima 
which  indicates  similar  trends. 

Of  greater  interest  is  the  case  for  a  cylinder  with  0  =  2  as  shown  in  Figure  7.6. 
The  Strouhal  number  increases  steadily  with  increasing  Reynolds  number  up  to  the 
point  where  Re  ~  600.  The  Strouhal  number  then  falls  sharply  and  maintains  a 
fairly  constant  value  over  the  remainder  of  the  Reynolds  number  range.  This  result 
is  in  good  agreement  with  experimental  results  which  indicated  the  same  behavior. 
The  prediction  of  the  point  where  the  discontinuity  occurred  was  also  very  close  to 
the  experimental  measurements. 

The  case  of  the  cylinder  with  0  =  3  shown  in  Figure  7.7  was  similar  to  the  case 
of  0  =  2.  Again,  the  Strouhal  number  increased  with  increasing  Reynolds  number 
up  to  a  value  of  Re  ^  1600.  At  this  point,  the  plot  then  changed  discontinuous!}' 
and  the  Reynolds  number  remained  reasonably  constant  over  the  remainder  of  the 
plot.  Qualitatively,  the  behavior  of  the  prediction  agrees  well  with  the  experimental 
data  provided  by  Okajima.  However,  the  discontinuity  in  the  plot  occurred  at  a 
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significantly  greater  Reynolds  number  than  cited  by  Okajima.  Nonetheless,  the 
prediction  of  the  behavior  of  the  plots  of  Strouhal  number  versus  Reynolds  number 
for  all  value  of  0  were  reasonably  accurate. 

Variation  in  Aerodynamic  Force  Coefficients 


A  variation  in  the  aerodynamic  force  coefficients  with  both  Reynolds  number 
and  aspect  ratio  was  indicated  by  the  data.  The  exhibited  trends  were  more  solidly 
established  for  drag,  with  the  variation  in  drag  with  both  Reynolds  number  and  0 
being  quite  evident. 

Figure  7.8  presents  the  variation  in  the  drag  coefficient  for  the  different  values 
of  0.  Clearly  the  drag  coefficient  decreases  for  increases  in  the  aspect  ratio.  This 
is  a  trend  which  was  indicated  by  the  study  of  Laneville  and  Yong.  A  variation  in 
drag  coefficient  with  the  Reynolds  number  was  also  apparent  as  shown  in  Figure 
7.9.  The  trend  is  most  clearly  seen  in  the  data  for  the  cylinder  with  0  =  2,  where 
the  drag  steadily  increases  with  increasing  Reynolds  number.  This  is  the  trend 
indicated  by  the  data  of  Davis  and  Moore.  However,  this  trend  was  more  subtle  for 
the  other  values  of  0. 

The  variations  in  the  lift  coefficient  with  0  shown  in  Figure  7.10  did  not  exhibit 
the  expected  behavior.  There  are  some  indications  that  the  lift  may  increase  slightly 
with  increasing  0,  but  the  trend  is  not  well  established.  Even  more  curious  was 
the  variation  in  the  lift  with  the  Reynolds  number  shown  in  Figure  7.11.  The  data 
in  this  case  clearly  indicate  that  the  lift  decreases  for  increases  in  the  Reynolds 
number.  Data  provided  by  Davis  and  Moore  indicates  a  relationship  between  the 
lift  and  Reynolds  number  inverse  to  that  derived  in  this  study. 
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FIGURE  7.4:  Absolute  Velocity  Vectors  for  /3 
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FIGURE  7.9:  Variation  in  the  Drag  Coefficient 
with  Reynolds  Number. 
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CHAPTER  VIII 
CONCLUSIONS 


This  study  proposed  to  demonstrate  the  ability  to  predict  the  phenomena  of 
vortex  shedding  from  various  two-dimensional  bodies.  For  every  case  tested  with 
the  exception  of  a  flat  plate  on  a  ground  plane,  periodic  vortex  shedding  was 
achieved. 

For  the  initial  case  involving  a  flat  plate  in  a  uniform  flow  (Chapter  V),  the 
results  were  favorable.  Predicted  Strouhal  numbers  agreed  well  with  experimental 
data.  A  calculation  of  pressures  in  the  wake  of  the  plate  indicated  a  region  behind 
the  plate  where  the  predicted  pressures  were  much  lower  than  the  experimental  data; 
however,  the  difference  may  be  due  to  difficulties  in  experimentally  measuring  the 
pressures.  Predictions  of  the  drag  also  indicated  deviations  from  the  drag  function 
provided  by  Roshko,  but  the  accuracy  of  RoshkoTs  function  was  unknown. 

With  respect  to  parameters  affecting  the  solution  of  the  finite-difference  equa¬ 
tions,  some  significant  results  were  discovered.  Refinement  of  the  finite-difference 
grid  yielded  only  a  slight  improvement  in  the  prediction  of  pressure  distributions, 
and  actually  hampered  the  prediction  of  Strouhal  numbers.  The  inclusion  of  turbu¬ 
lence  was  found  to  have  very  little  effect,  and  was  subsequently  excluded  from  the 


solution  process. 
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The  prediction  of  vortex  shedding  from  a  fiat  plate  near  a  ground  plane  also 
yielded  favorable  results.  The  comparison  of  predicted  Strouhal  numbers  with  an 
empirically  derived  function  demonstrated  the  ability  to  predict  general  trends  in 
the  shedding  frequency.  Quantitatively,  there  were  significant  differences  in  the 
prediction  and  the  empirical  function;  however,  this  may  be  attributable  to  the  fact 
that  the  prediction  was  two-dimensional  while  the  function  was  derived  for  a  three- 
dimensional  case.  Also,  the  expected  increase  in  the  aerodynamic  force  coefficients 
as  the  body  approached  the  ground  plane  were  accurately  predicted. 

Finally,  the  prediction  of  vortex  shedding  from  rectangular  cylinders  was  suc¬ 
cessful.  Strouhal  numbers  for  cylinders  of  differing  aspects  ratios  were  predicted 
with  reasonable  accuracy.  This  prediction  also  properly  included  the  discontinuity 
in  the  plots  of  Strouhal  number  versus  Reynolds  number  for  cylinders  with  aspect 
ratios  in  the  range  2  <  f3  <  3.  Unfortunately,  the  flow  pattern  that  gives  rise  to  this 
discontinuity  was  not  clearly  apparent.  Predictions  of  the  drag  coefficient  on  the 
cylinder  compare  well  with  previous  studies;  however,  the  lift  coefficient  behaves  in 
a  manner  contrary  to  prior  indications. 
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